require( googlesheets4 )
require(dplyr)
require(magrittr)
require(ggplot2)
require(reshape2)
require(cowplot)
if(dir.exists("output")!=TRUE) dir.create("output") # check if output directory out exists, if false createRepeability of Frog Morphometric Measurements
Assessing Repeatability of Morphometric Measurements
We want to measure repeatability of our morphometric measurements, both within individual measurers and between measurers. We plan to combine data with that measured by Julio Rivera in 2015, so we want to know that the data are compatible.
The accuracy we are aiming for is 95%, or 5% measurement error (both within and between individual measurers).
Each measurer measured the same 7 frogs (Hylophorbus sp. from Buyetai: JR306,308,311,320,321,324), which was repeated (2 sets of measurements per measurer).
Setup
Read in our data from a Google Spreadsheet, and convert measurer, session and jr number to factors.
file <- "https://docs.google.com/spreadsheets/d/1-w62GXvKwQ868dwiMVzgPdLYbgSkXxSqasE-YoaGOUc/edit#gid=0"
gs4_deauth() # not a private sheet, so no need for authentication
alldat <- as.data.frame(read_sheet(file)) %>% ### INPUT DATA from googlesheet
filter(!is.na(measurer)) # remove spacer rows✔ Reading from "Frog Repeatability Measurements".
✔ Range 'Sheet1'.
alldat %<>% mutate_at(c("measurer", "session", "jr"), as.factor) %>%
mutate( shape = case_when( measurer=="JR" ~ 25, measurer!="JR" ~ 19 )) %>%
mutate( shape = as.factor(shape))
jrblock <- alldat[alldat$measurer=="JR",] # copy JR data into all sessions (0,1,2)
jrblock2 <- jrblock0 <- jrblock
jrblock0$session <- 0
jrblock2$session <- 2
alldat <- rbind(alldat, jrblock0, jrblock2)Repeatability
This requires two full sets of measurements, so this is still in progress.
mod <- with(alldat, summary(aov(lm( svl ~ jr ))))
s2_within <- ms_within <- mod[[1]][2,3]
s2_within[1] 0.3160116
ms_among <- mod[[1]][1,3]
s2_among <- (ms_among-ms_within)/2
ME <- s2_within/(s2_within+s2_among) * 100
ME[1] 0.405364
Functions
These function creates the plots that are repeated for each morphometric variable.
dat <- alldat %>% filter(session==0)
p <- dat %>% ggplot(aes(svl, femur, color=measurer, label=jr))
v <- dat %>% ggplot(aes(jr, svl, group=jr, color=measurer, shape=measurer))
make_plots <- function(p, v) {
q1 <- p + geom_point(size = 3) +
geom_smooth( aes(group=measurer), method="lm", alpha=.1) +
geom_text(nudge_y = .15) +
theme_bw()
q2 <- p + geom_point(size = 3) +
geom_smooth( method="lm", alpha=.1) +
geom_text(nudge_y = 1) +
facet_grid( measurer ~ . ) +
theme_bw()
plot_grid(q1, q2, labels="AUTO")
}
make_violins <- function(v) {
v + geom_violin() +
geom_jitter(aes(x=jr, color=measurer, shape=measurer, size=measurer), width=.2) +
scale_shape_manual(values=c(19,19,17,19,19,19)) +
scale_size_manual(values=c(3,3,5,3,3,3))
}Set up the ggplots:
dat <- alldat %>% filter(session==0) # filter just for session 0 (naive)
femp <- dat %>% ggplot(aes(svl, femur, color=measurer, label=jr))
tibp <- dat %>% ggplot(aes(svl, tibiofibula, color=measurer, label=jr))
tarp <- dat %>% ggplot(aes(svl, tarsus, color=measurer, label=jr))
footp <- dat %>% ggplot(aes(svl, foot, color=measurer, label=jr))
hwp <- dat %>% ggplot(aes(svl, headW, color=measurer, label=jr))
hlp <- dat %>% ggplot(aes(svl, headL, color=measurer, label=jr))
hump <- dat %>% ggplot(aes(svl, humerus, color=measurer, label=jr))
radp <- dat %>% ggplot(aes(svl, radioulna, color=measurer, label=jr))
handp <- dat %>% ggplot(aes(svl, hand, color=measurer, label=jr))
sv <- dat %>% ggplot(aes(jr, svl, group=jr))
femv <- dat %>% ggplot(aes(jr, femur, group=jr))
tibv <- dat %>% ggplot(aes(jr, tibiofibula, group=jr))
tarv <- dat %>% ggplot(aes(jr, tarsus, group=jr))
footv <- dat %>% ggplot(aes(jr, foot, group=jr))
hwv <- dat %>% ggplot(aes(jr, headW, group=jr))
hlv <- dat %>% ggplot(aes(jr, headL, group=jr))
humv <- dat %>% ggplot(aes(jr, humerus, group=jr))
radv <- dat %>% ggplot(aes(jr, radioulna, group=jr))
handv <- dat %>% ggplot(aes(jr, hand, group=jr))Validating Landmarks (session 0)
Plot all of our data against the original dataset (JR in green), as well as faceted by individual. The goal here is to validate the landmarks.
SVL
make_violins(sv)Looks like everyoneʻs SVL are a little short. Basically JR304 and JR306 match well, but the others, JR308, JR311, JR320, JR321, JR322, JR324 are a little short. Did he press them down? Or use a ruler?
Diana and Ke investigated and found that gently pressing on the frog was necessary to reapeat JRʻs measurements. Just enough to maximize the length.
Femur
make_violins(femv)make_plots( femp )Looks like my femur measurements (MB) are a little short - JR must have measured from the vent rather than from the midline perpendicular to the femur.
The landmarks are from the edge of the knee to the midline at the vent.
Tibiofibula
make_violins(tibv)make_plots( tibp )Tibiofibula looks good, except for what looks like a typo in DGʻs data?
We measured the bone, from end to end.
Tarsus
make_violins(tarv)make_plots( tarp )Not sure what I did! Lol. Good example of the unfortunate outlier being the largest throwing the whole trend off. A few measurers are off on some specimens.
Foot
make_violins(footv)make_plots( footp )Looks good!
Head Width
make_violins(hwv)make_plots( hwp )Looks good! Measured at the widest part of the head or at the center of the tympanum (no further back than the tympanum).
Head Length
make_violins(hlv)make_plots( hlp )Oh oh. We are not using the right Head Length landmark.
Diana and Ke discovered the the measurement is from the tip of the snout to the center of the tympanum.
Humerus
make_violins(humv)make_plots( hump )Hmm. We are all a little consistently lower than JR, some more so.
The landmark for the end of the humerus can be most easily found by gently pressing the caliper onto the chest, which will move the arm. It is from the articulation of the humerus with the scapula to the edge of the elbow.
Radioulna
make_violins(radv)make_plots( radp)Need to work on this one too.
Hand
make_violins(handv)make_plots( handp )Looks good!
Conclusions
We need to confirm the landmarks JR used for the limb segment and head length measurements (foot and hand are OK). It looks promising for repeatability but we need to confirm after we get two full sets of measurements.
Session 1
Set up the ggplots:
dat <- alldat %>% filter(session==1)
femp <- dat %>% ggplot(aes(svl, femur, color=measurer, label=jr))
tibp <- dat %>% ggplot(aes(svl, tibiofibula, color=measurer, label=jr))
tarp <- dat %>% ggplot(aes(svl, tarsus, color=measurer, label=jr))
footp <- dat %>% ggplot(aes(svl, foot, color=measurer, label=jr))
hwp <- dat %>% ggplot(aes(svl, headW, color=measurer, label=jr))
hlp <- dat %>% ggplot(aes(svl, headL, color=measurer, label=jr))
hump <- dat %>% ggplot(aes(svl, humerus, color=measurer, label=jr))
radp <- dat %>% ggplot(aes(svl, radioulna, color=measurer, label=jr))
handp <- dat %>% ggplot(aes(svl, hand, color=measurer, label=jr))
sv <- dat %>% ggplot(aes(jr, svl, group=jr))
femv <- dat %>% ggplot(aes(jr, femur, group=jr))
tibv <- dat %>% ggplot(aes(jr, tibiofibula, group=jr))
tarv <- dat %>% ggplot(aes(jr, tarsus, group=jr))
footv <- dat %>% ggplot(aes(jr, foot, group=jr))
hwv <- dat %>% ggplot(aes(jr, headW, group=jr))
hlv <- dat %>% ggplot(aes(jr, headL, group=jr))
humv <- dat %>% ggplot(aes(jr, humerus, group=jr))
radv <- dat %>% ggplot(aes(jr, radioulna, group=jr))
handv <- dat %>% ggplot(aes(jr, hand, group=jr))SVL
make_violins(sv)Much better!
Femur
make_violins(femv)make_plots( femp )Better!
Tibiofibula
make_violins(tibv)make_plots( tibp )Better!
Tarsus
make_violins(tarv)make_plots( tarp )Better but Diana still a bit low and Ke a bit high
Foot
make_violins(footv)make_plots( footp )Still good!
Head Width
make_violins(hwv)make_plots( hwp )Good. I may have been too generous.
Head Length
make_violins(hlv)make_plots( hlp )Much better!
Humerus
make_violins(humv)make_plots( hump )Several folks are still a bit low.
Radioulna
make_violins(radv)make_plots( radp)Now looks like random error, still a bit high though.
Hand
make_violins(handv)make_plots( handp )Now looks like random error.
Session 2
dat <- alldat %>% filter(session==2)
femp <- dat %>% ggplot(aes(svl, femur, color=measurer, label=jr))
tibp <- dat %>% ggplot(aes(svl, tibiofibula, color=measurer, label=jr))
tarp <- dat %>% ggplot(aes(svl, tarsus, color=measurer, label=jr))
footp <- dat %>% ggplot(aes(svl, foot, color=measurer, label=jr))
hwp <- dat %>% ggplot(aes(svl, headW, color=measurer, label=jr))
hlp <- dat %>% ggplot(aes(svl, headL, color=measurer, label=jr))
hump <- dat %>% ggplot(aes(svl, humerus, color=measurer, label=jr))
radp <- dat %>% ggplot(aes(svl, radioulna, color=measurer, label=jr))
handp <- dat %>% ggplot(aes(svl, hand, color=measurer, label=jr))
sv <- dat %>% ggplot(aes(jr, svl, group=jr))
femv <- dat %>% ggplot(aes(jr, femur, group=jr))
tibv <- dat %>% ggplot(aes(jr, tibiofibula, group=jr))
tarv <- dat %>% ggplot(aes(jr, tarsus, group=jr))
footv <- dat %>% ggplot(aes(jr, foot, group=jr))
hwv <- dat %>% ggplot(aes(jr, headW, group=jr))
hlv <- dat %>% ggplot(aes(jr, headL, group=jr))
humv <- dat %>% ggplot(aes(jr, humerus, group=jr))
radv <- dat %>% ggplot(aes(jr, radioulna, group=jr))
handv <- dat %>% ggplot(aes(jr, hand, group=jr))SVL
make_violins(sv)Much better!
Femur
make_violins(femv)make_plots( femp )Better!
Tibiofibula
make_violins(tibv)make_plots( tibp )Better!
Tarsus
make_violins(tarv)make_plots( tarp )Better but Diana still a bit low and Ke a bit high
Foot
make_violins(footv)make_plots( footp )Still good!
Head Width
make_violins(hwv)make_plots( hwp )Good. I may have been too generous.
Head Length
make_violins(hlv)make_plots( hlp )Much better!
Humerus
make_violins(humv)make_plots( hump )Several folks are still a bit low.
Radioulna
make_violins(radv)make_plots( radp )Now looks like random error, still a bit high though.
Hand
make_violins(handv)make_plots( handp )All together
dat <- alldat
sv <- dat %>% ggplot(aes(jr, svl, group=jr))
femv <- dat %>% ggplot(aes(jr, femur, group=jr))
tibv <- dat %>% ggplot(aes(jr, tibiofibula, group=jr))
tarv <- dat %>% ggplot(aes(jr, tarsus, group=jr))
footv <- dat %>% ggplot(aes(jr, foot, group=jr))
hwv <- dat %>% ggplot(aes(jr, headW, group=jr))
hlv <- dat %>% ggplot(aes(jr, headL, group=jr))
humv <- dat %>% ggplot(aes(jr, humerus, group=jr))
radv <- dat %>% ggplot(aes(jr, radioulna, group=jr))
handv <- dat %>% ggplot(aes(jr, hand, group=jr))
wrap_violin <- function(v) {
v + geom_violin() +
geom_jitter(aes(x=jr, color=measurer, shape=measurer, size=measurer), width=.2) +
scale_shape_manual(values=c(19,19,17,19,19,19)) +
scale_size_manual(values=c(3,3,5,3,3,3)) +
facet_wrap(. ~ session )
}SVL
wrap_violin(sv)Femur
wrap_violin(femv)Tibiofibula
wrap_violin(tibv)Tarsus
wrap_violin(tarv)Foot
wrap_violin(footv)Head Width
wrap_violin(hwv)Head Length
wrap_violin(hlv)Humerus
wrap_violin(humv)Radioulna
wrap_violin(radv)Hand
wrap_violin(handv)